# robustness tests testing different thresholds for
# the agricultural variable

library(specr)
library(fixest)
library(readr)
library(broom.mixed)
library(texreg)
library(dplyr)
library(fixest)
library(patchwork)

# Cross sectional linear FE models that use aggregate data (2000-2019)

source("funs_and_constants.R")

# load in data
epr_yearly_group_clusters_h1 <- read_csv("data/epro_data_h1.csv")

mean_h1 <- epr_yearly_group_clusters_h1 %>% 
  filter(year >= 2000) %>%# only obs after 2000 because of EPR 2.0 coding
  group_by(gwgroupid, countries_gwid) %>%
  summarize(
    across(c(share_disagreement, disagreement, n_orgs, resource_and_agriculture_15, resource_and_agriculture_20,
             regaut, groupsize, incidence_flag, status_pwrrank, any_multiethnic, nightlight_total_pc_log, n_tek_groups),
           ~ mean(.x, na.rm = T)),
    none_one_or_both_nats_agri_char_15 = Mode(as.factor(none_one_or_both_nats_agri_char_15), na.rm = TRUE),
    none_one_or_both_nats_agri_char_20 = Mode(as.factor(none_one_or_both_nats_agri_char_20), na.rm = TRUE),
  )

full_model_h1_15 <- feols(disagreement ~ resource_and_agriculture_15 + n_orgs + regaut + groupsize + incidence_flag + status_pwrrank + any_multiethnic + nightlight_total_pc_log + n_tek_groups | countries_gwid,
                       cluster = ~countries_gwid, data = mean_h1 %>% as.data.frame())

full_model_h1_n_resources_15 <- feols(disagreement ~ none_one_or_both_nats_agri_char_15 + n_orgs + regaut + groupsize + incidence_flag + status_pwrrank + any_multiethnic + nightlight_total_pc_log + n_tek_groups | countries_gwid,
                                   cluster = ~countries_gwid, data = mean_h1)


full_model_h1_20 <- feols(disagreement ~ resource_and_agriculture_20 + n_orgs + regaut + groupsize + incidence_flag + status_pwrrank + any_multiethnic + nightlight_total_pc_log + n_tek_groups | countries_gwid,
                       cluster = ~countries_gwid, data = mean_h1 %>% as.data.frame())

full_model_h1_n_resources_20 <- feols(disagreement ~ none_one_or_both_nats_agri_char_20 + n_orgs + regaut + groupsize + incidence_flag + status_pwrrank + any_multiethnic + nightlight_total_pc_log + n_tek_groups | countries_gwid,
                                   cluster = ~countries_gwid, data = mean_h1)

screenreg(list(full_model_h1_15, full_model_h1_n_resources_15, full_model_h1_20, full_model_h1_n_resources_20),
          stars = stars,
          custom.coef.map = coef_map
          )

texreg(list(full_model_h1_15, full_model_h1_n_resources_15, full_model_h1_20, full_model_h1_n_resources_20),
       stars = stars,
       custom.coef.map = coef_map,
       custom.header = list("Threshold: 15 perc" = 1:2, "Threshold: 20 perc" = 3:4),
       table = FALSE,
       include.proj.stats = FALSE,
       file = "tables/robustness_agri_thresholds.tex")
